set more off

cd "Replication"

********************************************************************************
use "Data\iac_data.dta", clear
drop roads

merge m:1 dist61_id using "Data\census61.dta", nogen keep(match)

merge 1:1 dist61_id year using "Data\icrisat_data.dta", nogen keep(master matched) keepusing(rice_tai wheat_tai ptmrkt roads)

*Data errors
replace ywheat = . if ywheat == 15.02648 
replace yrice= . if yrice == 24

*Outcomes
gen ahyvwheat_gca = 100*ahyvwheat/gca
gen ahyvrice_gca = 100*ahyvrice/gca
gen prop_wheat = 100*ahyvwheat/awheat
gen prop_rice = 100*ahyvrice/arice
gen l_yrice = ln(yrice)
gen l_ywheat = ln(ywheat)
foreach var of varlist qnitro qphosphate qpotash {
gen `var'_gca = `var'/gca
gen `var'_nca = `var'/nca
}
gen rice_tai_area = 100*rice_tai/arice
gen wheat_tai_area = 100*wheat_tai/awheat
gen l_roads = ln(roads)
gen l_markets = ln(ptmrkt)

*Data errors
foreach var of varlist wheat_tai_area rice_tai_area prop_wheat prop_rice {
replace `var'=100 if `var' > 100 & `var' !=.
}

*****************************************************************************
*Covariates
gen lit_rate = 100*(m_lit+f_lit)/(tot_m_pop+tot_f_pop)
gen prop_rural = 100*t_r_pop/tot_pop

*Generating pre-GR values
foreach var of varlist arice awheat gca {
bys dist61: egen x = mean(`var') if year <=1966
bys dist61: egen `var'_pre = max(x)
drop x
}

gen arice_gca_pre = 100*arice_pre/gca_pre
gen awheat_gca_pre = 100*awheat_pre/gca_pre

gen x = wage if year == 1965
bys dist61_id: egen wage65 = min(x)
gen l_wage65 = ln(wage65)
drop x

********************************************************************************
merge m:1 dist61_id using "Data\radio_cov_data.dta", keep(matched)  nogen
merge 1:1 dist61_id year using "Data\radio_cov_share.dta", keep(master matched) nogen

sort dist61_id year
keep if year>=1966 & year <=1978

gen y = year-1966

foreach var of varlist awheat_gca_pre arice_gca_pre l_wage65 lit_rate prop_rural tot_pop_sq {
gen `var'_y1 = `var'*y
gen `var'_y2 = `var'*y^2
gen `var'_y3 = `var'*y^3
gen `var'_y4 = `var'*y^4
}

*State-year groups
drop state
encode statename, gen(state)
tab state, gen(state_)
egen state_year = group(state year)

foreach var of varlist state_1-state_12 {
gen `var'_y = `var'*y
gen `var'_y2 = `var'*y^2
gen `var'_y3 = `var'*y^3
gen `var'_y4 = `var'*y^4
}

*Radio coverage
gen reg_radio = (radial_radio_date_min<mdy(1,1,year))
replace reg_radio_share = 0 if reg_radio_share==.
gen vb_radio = vb_year<year

*Polynomial trend
global trend_poly "arice_gca_pre_y? awheat_gca_pre_y? state_*_y? lit_rate_y? prop_rural_y? tot_pop_sq_y? l_wage65_y?" 

*Defining panel data
xtset dist61_id year
sort dist61_id year

****Tables******
*Table 1
estpost tabstat gca nca ahyvwheat_gca ahyvrice_gca prop_wheat prop_rice wheat_tai_area rice_tai_area qnitro_nca qphosphate_nca qpotash_nca ywheat yrice roads ptmrkt if year == 1966, statistics(mean sd count) columns(statistics)
eststo m1
estpost tabstat gca nca ahyvwheat_gca ahyvrice_gca prop_wheat prop_rice wheat_tai_area rice_tai_area qnitro_nca qphosphate_nca qpotash_nca ywheat yrice roads ptmrkt if year == 1978, statistics(mean sd count) columns(statistics)
eststo m2
esttab m1 m2 using "Output/Table1.csv", replace cells("mean(fmt(3)) sd(fmt(3)) count(fmt(0))") nomtitle noobs nonumbers 

*Table 2
estpost tabstat reg_radio reg_radio_share vb_radio awheat_gca_pre arice_gca_pre wage65 lit_rate prop_rural tot_pop_sq if year == 1966, statistics(mean sd count) columns(statistics)
eststo m1
estpost tabstat reg_radio reg_radio_share vb_radio if year == 1978, statistics(mean sd count) columns(statistics)
eststo m2
esttab m1 m2 using "Output/Table2.csv", replace cells("mean(fmt(3)) sd(fmt(3)) count(fmt(0))") nomtitle noobs nonumbers 

*Table 3
eststo clear
foreach var of varlist ahyvwheat_gca prop_wheat ahyvrice_gca prop_rice {
eststo: xi: xtreg `var' reg_radio i.year $trend_poly , fe cluster(dist61)
estadd ysumm
eststo: xi: xtreg `var' reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
}
esttab using "Output/Table3.csv", replace keep(reg_radio reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N, fmt(3 0) labels("Adj. R^2" "Obs.")) 

*Table 4
eststo clear
foreach var of varlist ahyvwheat_gca prop_wheat ahyvrice_gca prop_rice  {
eststo: xi: xtreg `var' reg_radio F.reg_radio i.year $trend_poly , fe cluster(dist61)
estadd ysumm
test reg_radio = F.reg_radio
estadd scalar pval =  r(p)
eststo: xi: xtreg `var' reg_radio_share F.reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
test reg_radio_share = F.reg_radio_share
estadd scalar pval =  r(p)
}
esttab using "Output/Table4.csv", replace keep(reg_radio reg_radio_share F.reg_radio F.reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N pval, fmt(3 0 3) labels("Adj. $ R^2$" "Obs.")) 

*Table 5
eststo clear
foreach var of varlist ahyvwheat_gca prop_wheat ahyvrice_gca prop_rice {
eststo: xi: xtreg `var' vb_radio reg_radio i.year $trend_poly , fe cluster(dist61)
estadd ysumm
test reg_radio = vb_radio
estadd scalar pval =  r(p)
eststo: xi: xtreg `var' vb_radio reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
test reg_radio_share = vb_radio
estadd scalar pval =  r(p)
}
esttab using "Output/Table5.csv", replace keep(vb_radio reg_radio reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N pval, fmt(3 0 3) labels("Adj. R^2" "Obs.")) 

*Table 6
eststo clear
foreach var of varlist wheat_tai_area rice_tai_area {
eststo: xi: xtreg `var' reg_radio i.year $trend_poly , fe cluster(dist61)
estadd ysumm
eststo: xi: xtreg `var' reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
}
esttab using "Output/Table6.csv", replace keep(reg_radio reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N, fmt(3 0) labels("Adj. R^2" "Obs." ))

*Table 7
eststo clear
foreach var of varlist qnitro_gca qphosphate_gca qpotash_gca {
eststo: xi: xtreg `var' reg_radio i.year $trend_poly , fe cluster(dist61)
eststo: xi: xtreg `var' reg_radio_share i.year $trend_poly , fe cluster(dist61)
}
esttab using "Output/Table7.csv", replace keep(reg_radio reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N, fmt(3 0) labels("Adj. R^2" "Obs."))

*Table 8
eststo clear
foreach var of varlist l_roads l_markets {
eststo: xi: xtreg `var' reg_radio i.year $trend_poly , fe cluster(dist61)
estadd ysumm
eststo: xi: xtreg `var' reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
}
esttab using "Output/Table8.csv", replace keep(reg_radio reg_radio_share) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N, fmt(3 0) labels("Adj. R^2" "Obs.")) 

*Table 9
eststo clear
foreach var of varlist l_ywheat l_yrice {
eststo: xi: xtreg `var' reg_radio  i.year $trend_poly , fe cluster(dist61)
estadd ysumm
eststo: xi: xtreg `var' reg_radio_share i.year $trend_poly , fe cluster(dist61)
estadd ysumm
}
esttab using "Output/Table9.csv", replace keep(reg_radio reg_radio_share ) star(* 0.10 ** 0.05 *** 0.01) ///
wrap cells(b(fmt(3) star) se(par fmt(3))) collabel(none) label stats(r2_a N , fmt(3 0 ) labels("Adj. R^2" "Obs."))

****Figures****
gen reg_radio_year = yofd(radial_radio_date_min)
gen radio_year_ag = reg_radio_year
replace radio_year_ag = . if reg_radio_year>1978
replace radio_year_ag = 1965 if reg_radio_year<1966
gen y_ag_radio = year-radio_year_ag
replace y_ag_radio = -16 if radio_year_ag==.
recode y_ag_radio (-15/-8=-8) (8/13=8)
tab y_ag_radio, gen(y_ag_placebo)
recode y_ag_radio (-16 -8/-1 =0)
tab y_ag_radio, gen(y_ag)


*Figure 1
preserve
eststo clear
qui: eststo: xi: xtreg ahyvrice_gca y_ag_placebo3-y_ag_placebo18 i.year $trend_poly , fe cluster(dist61)
test y_ag_placebo3 = y_ag_placebo4 = y_ag_placebo5 = y_ag_placebo6 = y_ag_placebo7 = y_ag_placebo8  = y_ag_placebo9 = y_ag_placebo10
esttab using "Output/ahyvrice_year_wise_placebo.csv", replace keep(y_ag_placebo*) cells("b(fmt(3)) ci(par fmt(3))") 
insheet using "Output/ahyvrice_year_wise_placebo.csv", comma clear
drop if _n<=3
strparse v3, gen(ci) parse(",")
drop v1 v3
destring v2 ci1 ci2, ignore("[" "]") replace
drop if _n==17
egen year = seq()
replace year = year-8
rename v2 coeff
destring coeff, replace
twoway scatter coeff year, title(Rice) yline(0, lcolor(black)) connect(l) mcolor(black) legend(off) xlabel(-7/7 8 ">7", grid labsize(small)) ylabel(-4(0.5)7, angle(horizontal) labsize(small)) ysc(titlegap(1)) xsc(titlegap(2))  xtitle(Years after radio coverage, size(small)) ytitle(Effect on HYV Rice Area/GCA, size(small))|| scatter ci2 year, mcolor(black) connect(l) lcolor(black) lpattern(dash) || scatter ci1 year, mcolor(black) connect(l) lcolor(black) lpattern(dash)
graph export "Output/Figure5_1.png", replace
restore 
erase "Output/ahyvrice_year_wise_placebo.csv"

preserve
eststo clear
qui: eststo: xi: xtreg ahyvwheat_gca y_ag_placebo3-y_ag_placebo18 i.year $trend_poly , fe cluster(dist61)
test y_ag_placebo3 = y_ag_placebo4 = y_ag_placebo5 = y_ag_placebo6 = y_ag_placebo7 = y_ag_placebo8  = y_ag_placebo9 = y_ag_placebo10
esttab using "Output/ahyvwheat_year_wise_placebo.csv", replace keep(y_ag_placebo*) cells("b(fmt(3)) ci(par fmt(3))") 
insheet using "Output/ahyvwheat_year_wise_placebo.csv", comma clear
drop if _n<=3
strparse v3, gen(ci) parse(",")
drop v1 v3
destring v2 ci1 ci2, ignore("[" "]") replace
drop if _n==17
egen year = seq()
replace year = year-8
rename v2 coeff
destring coeff, replace
twoway scatter coeff year, title(Wheat) yline(0, lcolor(black)) connect(l) mcolor(black) legend(off) xlabel(-7/7 8 ">7", grid labsize(small)) ylabel(-4(0.5)7, angle(horizontal) labsize(small)) ysc(titlegap(1)) xsc(titlegap(2))  xtitle(Years after radio coverage, size(small)) ytitle(Effect on HYV Wheat Area/GCA, size(small))|| scatter ci2 year, mcolor(black) connect(l) lcolor(black) lpattern(dash) || scatter ci1 year, mcolor(black) connect(l) lcolor(black) lpattern(dash)
graph export "Output/Figure5_2.png", replace
restore 
erase "Output/ahyvwheat_year_wise_placebo.csv"

*Figure 6
preserve
eststo clear
eststo: xi: xtreg ahyvrice_gca y_ag2-y_ag9 i.year $trend_poly , fe cluster(dist61)
esttab using "Output/ahyvrice_year_wise.csv", replace keep(y_ag*) cells("b(fmt(3)) ci(par fmt(3))") 
insheet using "Output/ahyvrice_year_wise.csv", comma clear
drop if _n<=3
strparse v3, gen(ci) parse(",")
drop v1 v3
destring v2 ci1 ci2, ignore("[" "]") replace
drop if _n==9
egen year = seq()
rename v2 coeff
destring coeff, replace
twoway scatter coeff year, title(Rice) yline(0, lcolor(black)) connect(l) mcolor(black) legend(off) xlabel(1/7 8 ">7", grid labsize(small)) ylabel(-1(0.5)6, angle(horizontal) labsize(small)) ysc(titlegap(1)) xsc(titlegap(2))  xtitle(Years after radio coverage, size(small)) ytitle(Effect on HYV Rice Area/GCA, size(small))|| scatter ci2 year, mcolor(black) connect(l) lcolor(black) lpattern(dash) || scatter ci1 year, mcolor(black) connect(l) lcolor(black) lpattern(dash)
graph export "Output/Figure6_1.png", replace
restore 
erase "Output/ahyvrice_year_wise.csv"

preserve
eststo clear
eststo: xi: xtreg ahyvwheat_gca y_ag2-y_ag9 i.year $trend_poly , fe cluster(dist61)
esttab using "Output/ahyvwheat_year_wise.csv", replace keep(y_ag*) cells("b(fmt(3)) ci(par fmt(3))") 
insheet using "Output/ahyvwheat_year_wise.csv", comma clear
drop if _n<=3
strparse v3, gen(ci) parse(",")
drop v1 v3
destring v2 ci1 ci2, ignore("[" "]") replace
drop if _n==9
egen year = seq()
rename v2 coeff
destring coeff, replace
twoway scatter coeff year, title(Wheat) yline(0, lcolor(black)) connect(l) mcolor(black) legend(off) xlabel(1/7 8 ">7", grid labsize(small)) ylabel(-1(0.5)6, angle(horizontal) labsize(small)) ysc(titlegap(1)) xsc(titlegap(2))  xtitle(Years after radio coverage, size(small)) ytitle(Effect on HYV Wheat Area/GCA, size(small))|| scatter ci2 year, mcolor(black) connect(l) lcolor(black) lpattern(dash) || scatter ci1 year, mcolor(black) connect(l) lcolor(black) lpattern(dash)
graph export "Output/Figure6_2.png", replace
restore 
erase "Output/ahyvwheat_year_wise.csv"



